/**********************************
**** This file will read in the public use data set for the
	2019 Survey of American Science teachers and reproduce tables
	published in:
	
		Branch, Glenn, Ann Reid and Eric Plutzer. 2021. 
		"Teaching Evolution in U.S. Public Middle Schools: 
		Results of the First National Survey." 
		Evolution: Education & Outreach 14 (8) . 
		https://doi.org/10.1186/s12052-021-00145-z 
*/


***********************************************************
**** NOTE: TABLE COMMANDS WILL NOT RUN IN STATA 17
****       To replicate published tables, the version command
****       must be executed.
****       This file will run without error if executed in its
****       entirety.  Otherwise, the version command must be 
****       executed before each table request.

version 16

*use "2019_Survey_of_American_Science_Teachers_PUS_dataset.dta"

use "C:\Users\exp12\OneDrive - The Pennsylvania State University\x New projects\2019 Teacher Survey\0 Methodology\Replication materials\2019_Survey_of_American_Science_Teachers_PUS_dataset.dta"



* Eliminate HS Chemistry/Other teachers and MS teachers not teaching science
drop if flag_combined == 0

* Declare data as survey data
svyset [pw=final_weight]

*********************************************************
***** TOTAL HOURS
gen total_hours = q1_num_gen + q1_num_humanevo
replace total_hours = q1_num_gen if q1_num_humanevo == .
replace total_hours = q1_num_humanevo  if  q1_num_gen == .

*********************************************************
*** MESSAGES EMPHASIZED IN CLASS - BUILD TYPOLOGY

* 1 generate ordinal measures
gen creat_ordinal_1 = 0
replace creat_ordinal_1 = 1 if q2_evo_intelligentdesign == 2 // agree
replace creat_ordinal_1 = 2 if q2_evo_intelligentdesign == 1 // agree strongly

gen creat_ordinal_2 = 0
replace creat_ordinal_2 = 1 if q2_evo_idvalidalt == 2 // agree
replace creat_ordinal_2 = 2 if q2_evo_idvalidalt == 1 // agree strongly

egen creat_ordinal = rowmax(creat_ordinal_1 creat_ordinal_2)

***  labeling
label define ord  0 "Agrees with neither statement"
label define ord  1 "Agrees with at least one", add
label define ord  2 "Strongly agrees with at least one", add
label values creat_ordinal ord

*
gen evo_ordinal_1 = 0
replace evo_ordinal_1 = 1 if q2_evo_consensus == 2 // agree
replace evo_ordinal_1 = 2 if q2_evo_consensus == 1 // agree strongly

gen evo_ordinal_2 = 0
replace evo_ordinal_2 = 1 if q2_evo_theme == 2 // agree
replace evo_ordinal_2 = 2 if q2_evo_theme == 1 // agree strongly

gen evo_ordinal = evo_ordinal_1 + evo_ordinal_2


tab evo_ordinal_1 creat_ordinal if level == 2 & hs_chem_other==0 [aw=final_weight], missing cell



gen evo_free_discussion = 0
replace evo_free_discussion = 1 if q3_discussion == 1
replace evo_free_discussion = . if q3_discussion == .


*** PROFILES
gen new_typology = .
replace new_typology = 1 if (evo_ordinal_1 ==1 | evo_ordinal_1 ==2)  & ///
	creat_ordinal ==0
replace new_typology = 2 if (evo_ordinal_1 ==1 | evo_ordinal_1 ==2) & ///
	(creat_ordinal ==1 | creat_ordinal ==2) 
replace new_typology = 3 if evo_ordinal_1 ==0 & ///
	creat_ordinal ==0
replace new_typology = 4 if evo_ordinal_1 ==0 & ///
	(creat_ordinal ==1 | creat_ordinal ==2) 

*** Setup for labeling
label define prof  1 "Evolution as settled science"
label define prof  2 "Mixed messages", add
label define prof  3 "Avoidance", add
label define prof  4 "Creationism", add
label values new_typology prof


*** Flag for those sayin "not applicable" to all three items in the typology
* MS: 197 and HS: 57
gen all_na = 0
replace all_na = 1 if ///
	q2_evo_intelligentdesign == 5 & ///
	q2_evo_idvalidalt == 5 & ///
	q2_evo_consensus == 5


* Create flag for typology
gen flag_typology = 1
replace flag_typology = 0 if total_hours == 0 & q1_idcreation == 1 // Zero hours for both evo & creationism (N=248)
replace flag_typology = 0 if total_hours == . & q1_idcreation == . // Both missing (N=39)
replace flag_typology = 0 if total_hours == . & q1_idcreation == 1 // One missing and other zero (N=2)
replace flag_typology = 0 if total_hours == 0 & q1_idcreation == . // One missing and other zero (N=1)


************************************
** PERSONAL BELIEFS ABOUT EVOLUTION, CONSENSUS & BIBLE

gen young_earther = .
replace young_earther = 1 if q16_origin == 3
replace young_earther = 0 if q16_origin <= 2

gen evo_80_percent = .
replace evo_80_percent = 1 if q17_evol_consensus == 5
replace evo_80_percent = 0 if q17_evol_consensus <= 4

gen literal = .
replace literal = 0 if q31_bible == 2 | q31_bible == 3
replace literal = 1 if q31_bible == 1 


***********************************
*** Advanced science training
gen degree_in_science = 0
replace degree_in_science = 1 if q21_bachscience ==1 | q21_mastersci ==1 | q21_phdscience == 1

gen grad_in_science = 0
replace grad_in_science = 1 if q21_mastersci ==1 | q21_phdscience == 1

gen grad_education = 0
replace grad_education = 1 if q21_mastered ==1 | q21_doctored

*******************
* Seniority in three categories

gen seniority_x = .
replace seniority_x = 1 if q22 < 10
replace seniority_x = 2 if q22 >=10 & q22<20
replace seniority_x = 3 if q22 >=20 & q22<50

label define sn3 1 "Under 10 yrs seniority"
label define sn3 2 "10{c 150}19 yrs seniority", modify
label define sn3 3 "20+ yrs seniority", add
label values seniority_x sn3

* Pre-service preparation
****** COLLEGE CLASS IN EVOLUTION ******
svy: tabulate level q27_evo_college    ,  percent format(%9.2f) stubw(33) row 
svy: tabulate level q27_sectionevo_college    ,  percent format(%9.2f) stubw(33) row

replace q27_evo_college = 1 if q27_evo_college == .
replace q27_sectionevo_college = 1 if q27_sectionevo_college == .
replace q27_sectionevo_college = 1 if q27_sectionevo_college == .
replace q27_evo_conted = 1 if q27_evo_conted == .

* Improving labels
label variable degree_in_science   "Holds degree in science"	
label variable q27_sectionevo_colleg   "College courses with some evolution"	
label variable q27_evo_college   "College courses focused on evolution"	

* Combining course variables
gen course_prep = 0
replace course_prep = 1 if (q27_sectionevo_college -1) == 1 | (q27_evo_college-1) == 1
replace course_prep = 2 if ((q27_sectionevo_college -1)) + ((q27_evo_college-1)) >= 2
replace course_prep = 3 if ((q27_sectionevo_college -1)) + ((q27_evo_college-1)) >= 3

label define cp 0 "No courses in evolution", replace
label define cp 1 "One course covering evolution", add
label define cp 2 "Two courses covering evolution", add
label define cp 3 "Three or more courses", add
label values course_prep cp
label variable course_prep   "College coursework in evolution"


* NGSS hours

table ngss_level   if level==1 [aweight = final_weight], contents(mean total_hour freq ) row format(%9.2f) stubwidth(30)

* check significance
svy: reg total_hour	 i.ngss_level if level==1
reg total_hour	 i.ngss_level if level==1

***********************************************
***  TABLES exported to Excel with tabout  ****
***  Requires installation of TABOUT version 3
***  Available from tabout.net.au 
***  
***********************************************

* Designate directory to write out tables in Excel format
cd "...tables"

* Table 1 - Preparation

tabout course_prep level [aw=final_weight] ///C:\Users\exp12\OneDrive - The Pennsylvania State University\x New projects\
using paper2_T1_preservice.xlsx,  ///
replace location (3 1)style(xlsx) font(plain)  noborder ///
lwidth (13) cwidth (10) indent (1) ///
c(col)  f(1 ) per  seb ///
npos(col) nlab(N) ///
title(Number of evo courses )

svy: tabulate course_prep level ,  count stubwidth(30) format(%9.0f)

*** Detailed tables for appendix

* Collapse 3 and 4 or more for reporting:
clonevar evo_college_3 = q27_evo_college
replace  evo_college_3 = 4 if  evo_college_3 == 5
tab  evo_college_3 q27_evo_college, nolabel

clonevar evo_section_3 = q27_sectionevo_college
replace evo_section_3 = 4 if evo_section_3 == 5
tab evo_section_3 q27_sectionevo_college, nolabel

clonevar evo_conted_3 = q27_evo_conted
replace evo_conted_3 = 4 if evo_conted_3 == 5
tab evo_conted_3 q27_evo_conted, nolabel

tabout level evo_college_3 [aw=final_weight] ///
using paper2_TA7_preservice.xlsx,  ///
replace location (3 1)style(xlsx) font(plain)  noborder ///
lwidth (13) cwidth (10) indent (1) ///
c(row)  f(1 ) per  seb ///
npos(col) nlab(N) ///
title(Number of evo courses )

tabout level evo_section_3 [aw=final_weight] ///
using paper2_TA7_preservice.xlsx,  ///
append location (12 1)style(xlsx) font(plain)  noborder ///
lwidth (13) cwidth (10) indent (1) ///
c(row)  f(1 ) per  seb ///
npos(col) nlab(N) 

tabout level evo_conted_3 [aw=final_weight] ///
using paper2_TA7_preservice.xlsx,  ///
append location (19 1)style(xlsx) font(plain)  noborder ///
lwidth (13) cwidth (10) indent (1) ///
c(row)  f(1 ) per  seb ///
npos(col) nlab(N) 

* Table 2 - Consensus

tabout q17_evol_consensus level  [aw=final_weight] ///
using paper2_T2_consensus.xlsx,  ///
replace location (3 1)style(xlsx) font(plain)  noborder ///
lwidth (13) cwidth (10) indent (1) ///
c(col)  f(1 ) per  seb ///
npos(col) nlab(N) ///
title(Perceived consensus among scientists)

* Table 3 - Personal views

tabout q16 level  [aw=final_weight] ///
using paper2_T3_personal.xlsx,  ///
replace location (3 1)style(xlsx) font(plain)  noborder ///
lwidth (13) cwidth (17) indent (1) ///
c(col)  f(1 ) per  seb ///
npos(col) nlab(N) ///
title(Personal opinions on evolution)

* Table 4 - hours

tabout level q1_generalevo [aw=final_weight] ///
using paper2_table4.xlsx,  ///
replace location (3 1) style(xlsx) font(plain)  noborder ///
lwidth (13) cwidth (7) indent (1) ///
c(row)  f(1 ) per  seb ///
npos(col) nlab(N) ///
title(Number of hours devoted to evolution)

tabout level q1_humanevo [aw=final_weight] ///
using paper2_table4.xlsx,  ///
append location (12 1) style(xlsx) font(plain) noborder  ///
lwidth (13) cwidth (7) indent (1) ///
c(row)  f(1 ) per  seb ///
npos(col) nlab(N) 

* Table 5 - evo as fact

tabout level q2_evo_consensus [aw=final_weight] ///
using paper2_table5.xlsx,  ///
replace location (3 1) style(xlsx) font(plain)  noborder ///
lwidth (13) cwidth (10) indent (1) ///
c(row)  f(1 ) per  seb ///
npos(col) nlab(N) ///
title(Teacher emphasis when teaching evolution)

tabout level q2_evo_consensus [aw=final_weight] ///
using paper2_table5.xlsx if total_hours > 0 & q2_evo_consensus<5,  ///
append location (12 1) style(xlsx) font(plain)  noborder ///
lwidth (13) cwidth (10) indent (1) ///
c(row)  f(1 ) per  seb ///
npos(col) nlab(N) 

* For footnote
gen related_topics = 0
replace related_topics = 1 if q1_any_cell == 1 | q1_any_biodiv == 1

tab related q2_evo_consensus [aw=final_weight] if level==1, col
tab related q2_evo_consensus [aw=final_weight] if level==2, col
tab related q2_evo_consensus [aw=final_weight] , col

svy: tab total_hours q2_evo_consensus  if level==1, col percent

* Table 6 - hours on creationism

tabout level q1_idcreation [aw=final_weight] ///
using paper2_table6.xlsx,  ///
replace location (3 1) style(xlsx) font(plain) noborder  ///
lwidth (13) cwidth (7) indent (1) ///
c(row)  f(1 ) per  seb ///
npos(col) nlab(N) ///
title(Number of hours devoted to creationism)

tabout level q1_idcreation [aw=final_weight] ///
using paper2_table6.xlsx if total_hours > 0,  ///
append location (12 1) style(xlsx) font(plain) noborder  ///
lwidth (13) cwidth (7) indent (1) ///
c(row)  f(1 ) per  seb ///
npos(col) nlab(N) ///
title(For those who taught one or more hours of evolution)

* Table 7 - creationism emphasis

tabout level q2_evo_intell [aw=final_weight] ///
using paper2_table7.xlsx,  ///
replace location (3 1) style(xlsx) font(plain)  noborder ///
lwidth (13) cwidth (10) indent (1) ///
c(row)  f(1 ) per  seb ///
npos(col) nlab(N) ///
title(I emphasize intelligent design is valid alternative)

tabout level q2_evo_intell if q2_evo_intell<5 [aw=final_weight] ///
using paper2_table7.xlsx  ,  ///
append location (12 1) style(xlsx) font(plain)  noborder ///
lwidth (13) cwidth (10) indent (1) ///
c(row)  f(1 ) per  seb ///
npos(col) nlab(N)

tabout level q2_evo_idvalidalt [aw=final_weight] ///
using paper2_table7.xlsx ,  ///
append location (21 1) style(xlsx) font(plain)  noborder ///
lwidth (13) cwidth (10) indent (1) ///
c(row)  f(1 ) per  seb ///
npos(col) nlab(N) ///
title(I emphasize scientists view intelligent design as valid)

tabout level q2_evo_idvalidalt if q2_evo_idvalidalt<5 [aw=final_weight] ///
using paper2_table7.xlsx ,  ///
append location (30 1)  style(xlsx) font(plain)  noborder ///
lwidth (13) cwidth (10) indent (1) ///
c(row)  f(1 ) per  seb ///
npos(col) nlab(N)



*** TABLE 8 - Typology by level
* 8A. Full sample
tabulate new_typology level // unweighted N
svy: tabulate new_typology level, count stubwidth(30) // weighted N
svy: tabulate new_typology level, col stubwidth(30) percent format(%9.1f) 

* 8B. Typology only for those who teach 1+ hours to either evolution or creationism
tabulate new_typology level if flag_typology == 1 // unweighted N
svy: tabulate new_typology level if flag_typology == 1 , count stubwidth(30) //weighted N
svy: tabulate new_typology level if flag_typology == 1, col stubwidth(30) percent format(%9.1f) 



*************************************************
************         FIGURES       **************

*** Setup for labeling
label define ng  0 "Non-Framework", replace
label define ng  1 "Framework", add
label define ng  2 "Adopted NGSS", add
label values ngss_level ng
tab ngss_level

label define dgs  0 "No degree in science"
label define dgs  1 "One or more degrees in science", add
label values degree_in_science dgs


************ FIGURE 1 *************
svy: reg total_hour	 i.ngss_level  ///
	i.course_prep i.degree_in_science i.seniority_x ///
	if level==1 , cformat(%9.3f)	
	
coefplot , drop(_cons) msymbol(O) base grid(none) ///
mcolor(blue) lcolor(blue) levels(90) scheme(s1mono) ///
order(*ngss_level . *course_prep . *degree_in_science . *seniority_x) ///
xline(0, lcolor(red%40)) /// 
xtitle("   " /// 
"Fewer                                                            More " "hours                                                            hours")

************ FIGURE 2 *************
* Binary comparison
gen settled = 0
replace settled = 1 if new_typology == 1


svy: logit settled	 i.ngss_level  ///
	i.course_prep i.degree_in_science i.seniority_x ///
	if level==1 & flag_typology==1, cformat(%9.3f) or nofvlabel
	
testparm i.ngss_level

coefplot , drop(_cons) msymbol(0) base grid(none) eform(*) ///
order(*ngss_level . *course_prep . *degree_in_science . *seniority_x ) ///
mcolor(blue) lcolor(blue) levels(90) scheme(s1mono) xline(1, lcolor(red%40)) ///
xtitle("   " "Relative risk ratio (Odds ratio)" )


************ APPENDIX FIGURE 3 *************

svy: mlogit new_typology  i.ngss_level  ///
	i.course_prep i.degree_in_science i.seniority_x  ///
	if level==1 & flag_typology == 1, cformat(%9.3f) rrr

coefplot ///
        (,keep(*Mi*:*) drop(_cons) label(Mixed Message versus Settled Science) msymbol(D) mcolor(gray) lcolor(gray)) ///
        (,keep(*Av*:*) drop(_cons) label(Avoidance versus Settled Science) msymbol(O) mcolor(blue) lcolor(blue)) ///
		(,keep(*Cre*:*) drop(_cons) label(Creationism versus Settled Science) msymbol(S) mcolor(black) lcolor(black)) ///
        , eform xline(1, lcolor(red%40)) legend(rows(3)) levels(90) scheme(s1mono) grid(none) ///
		xtitle("   " "Relative Risk Ratio (Odds ratio)" )

*************************


